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ABSTRACT 

A broad swathe of astrophysical phenomena, ranging from tubular planetary nebulae 
C^ I through Herbig-Haro objects, radio-galaxy and quasar emissions to gamma-ray bursts 

and perhaps high-energy cosmic rays, may be driven by magnetically-dominated jets 

emanating from accretion disks. We give a self-contained account of the analytic theory 

^ ■ of non-relativistic magnetically dominated jets wound up by a swirling disk and making 

^^ I a magnetic cavity in a background medium of any prescribed pressure, p{z). We solve 

Cn ■ the time-dependent problem for any specified distribution of magnetic flux P(i?,0) 

'^ I emerging from the disk at z = 0, with any specified disk angular velocity fldiR). 

The physics required to do this involves only the freezing of the lines of force to the 
\^ I conducting medium and the principle of minimum energy. 

f^ . In a constant pressure environment the magnetically dominated cavity is highly 

collimated and advances along the axis at a constant speed closely related to the 
rS . maximum circular velocity of the accretion disk. Even within the cavity the field is 

strongly concentrated toward the axis. The twist in the jet's field < B^ > / < \Bz\ > 
is close to v^ and the width of the jet decreases upwards. By contrast when the back- 
"^ ' ground pressure falls off with height with powers approaching z^'* the head of the jet 

C^ . accelerates strongly and the twist of the jet is much smaller. The width increases to 

give an almost conical magnetic cavity with apex at the source. Such a regime may 

be responsible for some of the longest strongly collimated jets. When the background 

^% ' pressure falls off faster than z~^ there are no quasi-static configurations of well twisted 

fields and the pressure confinement is replaced by a dynamic effective pressure or a 
relativistic expansion. In the regimes with rapid acceleration the outgoing and incom- 
ing fields linking the twist back to the source are almost anti-parallel so there is a 
possibility that magnetic reconnections may break up the jet into a series of magnetic 
'smoke-rings' travelling out along the axis. 
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1 INTRODUCTION 

1.1 Orders of magnitude — the voltages generated 

We consider accretion disks of bodies in formation when the differential rotation drags around magnetic field lines. Although 
the moving magnetic fields inevitably generate electric fields, the resulting EMFs are perpendicular to the magnetic fields 
in the perfect conduction approximation. Such EMFs do not accelerate particles to high energies. However the world is not 
a perfect place; in regions where perfect MHD predicts very high current densities there may be too few charge carriers to 
carry those currents. In those regions perfect conductivity is not a good approximation and the fields are modified to allow 
an electric field component along the magnetic field so that the larger current can be generated by the few charge carriers 
available. With this background idea it becomes interesting to get rough estimates of EMFs that are likely to be around 
whatever their direction. It turns out that these EMFs scale with (v/c) where v is the maximum velocity in the accretion 
disk and is independent of the size or mass of the system. Relativistic systems in formation, however small, can generate EMFs 
whose voltages at least match those of the highest energy cosmic rays. However the timescales over which these systems persist 
and the total energies available over those time scales do of course depend on the mass of the system. Very crude estimates of 
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the EMFs can be made as follows. Let us suppose that there is a central object, a star or a black hole, of mass M surrounded 
by an accretion disk of mass C,M. Further let us suppose that a small fraction r) of the binding energy of the accretion disk 
is converted into magnetic field energy. Letting R be the inner radius of the disk we put (47r/3)i?^(i3'^/(87r)) = GC^rjM^ / R so 
B = {GC^rj / GY'^ v^ / R., where v^ — GM/R. Fields of a few hundred Gauss are found in T Tauri stars, so putting R — lO^^cm 
V — 300km/s we find (^r] — \Q^ . We shall use this dimensionless number to make estimates elsewhere. Now the EMF 
generated around a circuit that goes up the axis and back to the disk around a field line and finally along the disk back to the 
centre is about {v/c)BR in esu and 300 times that in volts so EMF= 300 {&C,r) / Cf^ {v / cf c^ = («/c)^10^^ vohs. Kronberg 
et al (2004) considered such accelerative processes in radio galaxies. 

Again only a fraction of this voltage will be available to accelerate particles; in the exact fiux freezing case it is ALL 
perpendicular to the fields, nevertheless it is an estimate of what might be available where the perfect conductivity approxi- 
mation breaks down. Although the argument is crude the answer is interesting in that it suggests that cosmic rays generated 
in microquasars may reach the same individual particle energies as those generated in quasars and radio galaxies which may 
match the highest cosmic ray energies. 

1.2 Jets 

Curtis photographed the jet in M 87 from the Lick Observatory in 1918. While it was soon found to be blue, the emission 
process was not understood despite Schott's detailed calculations of synchrotron radiation in his 1912 book. Finally after 
Shklovskii introduced this emission mechanism into astrophysics Baade (1956) showed the jet to be highly polarized which 
clinched it. Although Ryle's Cambridge group (1968) found many double radio lobes around large galaxies the next obvious 
jet came with the identification of the first quasar 3C273 in which the dominant radio source is not at the nucleus but at the 
other end of the optical jet (Hazard, Mackey & Shimmins 1963; Schmidt 1963). There were difficulties in understanding the 
powering of the radio lobes of galaxies. If all the energy were present as the lobes expanded outwards there would be more 
bright small ones. Finally this led Rees (1971) to suggest that the lobes must be continuously powered by as-yet-unseen jets 
feeding energy into the visible lobes. As radio astronomy moved to higher frequencies with greater resolution these jets duly 
appeared in both radio galaxies and quasars. All the above jets have dynamically significant magnetic fields seen via their 
synchrotron radiation; a particularly fine study of one, Herculis A, is found in Gizani & Leahy (2003). Magnetic fields are 
less obvious in the Herbig-Haro objects first found in star-forming regions in 1951. However it took the development of good 
infra-red detectors before the heavy obscuration was penetrated in the 1980s to reveal the jets feeding the emission. These jets 
have velocities of one or two hundred km/s, far less than the 0.1c — c speeds of the extragalactic jets. The jets around young 
stars are seen to be perpendicular to the accretion discs that generate them. Since the giant black-hole accretion disk theory 
of quasars Salpeter (1964), Lynden-BeU (1969), Bardeen (1970), Lynden-Bell & Rees (1971) Shakura & Sunyaev (1973, 1976) 
many have come to believe that the radio galaxies and quasars likewise have jets perpendicular to their inner disks. However 
it was not until the wonderful work on megamasers (Miyoshi et al. 1995) that such inner disks around black holes were 
definitively confirmed. In 1969 when I predicted that they would inhabit the nuclei of most major galaxies including our own, 
M31, M32, M81, M87, etc., the idea was considered outlandish, but now most astronomers take it for granted. Fine work by 
Kormendy (1995) and others on external galaxies and the beautiful results of Genzel (2003) and Chez (2004) on our own has 
totally transformed the situation. Meanwhile many jets have been found in objects associated with dying stars, SS433 and the 
micro- quasars Mirabel & Rodriguez (1999) being prominent examples within the Galaxy. Less energetic but more beautiful 
examples may be the tubular planetary nebulae that are associated with accretion discs of central mass- exchanging binary 
stars. These were brought to my attention by Mark Morris and recent evidence indicates that magnetism is important here 
too (Vlemmings et al 2006). Much more spectacularly the 7-ray bursts are now thought to come from accretion-disks around 
black holes within some supernova explosions. Poynting flows of electromagnetic energy appear to be one of the best ways 
of extracting the energy from beneath the baryons that would otherwise absorb the 7-rays, Uzdensky & MacFadyen (2006). 
Jets and coUimated outflows have also been invoked for giant stars. Remarkable examples are R Aquarii (Michalitsanos et al. 
1988) and the Egg nebula (Cohen et al. 2004) and the coUimated outflow in IRG 10011 (Vinkovic 2004). 

The systematic features of these diverse objects are that an accretion disk is present and the jet emerges along the rotation 
axis of the inner disk. Magnetic fields are important in the radio objects and may be important in all. The jet velocities are 
strongly correlated with the escape velocities and therefore with the circular velocities in the disk close to the central object. 
In very collapsed objects these velocities are relativistic but in star-forming regions the jet velocities are less than lO^'^c. The 
obvious similarity of the jet structures and collimation despite such velocity differences strongly suggest that relativity is not 
a determining factor in the making of these jets. The thesis that I put forward in 1996 (paper II) and repeat here, is that 
all these jets are magnetically driven, the common feature being that the Poynting flux dominates the energy transport in 
the jet. However I do not exclude the possibility of material being entrained as the jet makes its way through the material 
that surrounds it. Nevertheless we consider that magnetism is the driver and the prime reason for collimation is the magnetic 
twisting combined with a weak external pressure which is dynamically enhanced by inertia. It is the electromagnetic field 
that generates the relativistic particles in radio jets so, even if their total energies evolve to become comparable, the magnetic 
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energy comes first and drives the whole phenomenon. While this paper is concerned with jets from systems with accretion 
disks, the jets near pulsars are probably closely related. There the magnetic field is rotated by the neutron star but the weaker 
field at large radii may be heavily loaded by the inertia close to the light cylinder. Field lines may be twisted from below far 
faster than they can rotate across the light cylinder. The resulting twisting up of the field within the cylinder may result in 
jet phenomena with similarities to those described here. 

The problem of jet coUimation was emphasised by Wheeler (1971) at the Vatican Conference on the Nuclei of Galaxies 
in 1970. He drew attention to the computations of Leblanc & Wilson (1970) who found a remarkable jet generated on the 
axis of a rotating star in collapse. This may be hailed as the first gamma-ray burst calculation and the magnetic cavities 
discussed here are in essence analytic calculations based on the same mechanism. We have concentrated on the simplest case 
of force-free magnetic fields within the cavity. In earlier work others produced good coUimation by twisting up a magnetic 
field that was imposed from infinity. An early paper on this was by Lovelace (1976) and the rather successful simulations 
of Shibata & Uchida (1985, 1986) are based on this theme. I consider that the imposition of a straight field from infinity 
ducks the question of why coUimation exists. The straight field imposes it. Much work has concentrated on the hard problem 
of winds carrying a significant mass fiux from a rotating star. This subject is well covered in Mestel's book (1999) Li et al 
(2001), Lovelace et al (2002) and Sakarai (1987) found a weak asymptotic coUimation and Heyvaerts & Norman (2002) have 
recently concluded a thorough study of the asymptotics of wind coUimation. Bogovalov & Tsinganos (2003) have tackled the 
difficulties of coUimating mass loaded fiows from central objects and give models with the near-axis part of the fiow well 
coUimated despite the centrifugal force. Blandford & Payne (1982) discuss winds launched by centrifugal force. Lovelace & 
Romanova (2003) Li et al (2006) made numerical calculations based upon the differential winding we proposed in Papers I 
& II. The stability problems of twisted jets have been tackled in both the linear and non-linear regimes by Appl, Lery & 
Baty (2000) and Lery, Baty & Appl (2000) while Thompson, Lyutikov and Kulkarni (2002) have applied to magnetars the 
self-similar fields found in Paper I. In section 5 we find that some magnetic cavities fioat upwards like bubbles thus fulfilling 
the ideas of Gull & Northover (1976). The laboratory experiments of Lebedev et al. (2005) give some support for the type of 
models given here. A fine review of extragalactic jets was given by Begelman et al. (1984). See also the more recent work of 
Pudritz et al 2006, Ouyed et al 2003. 

1.3 Outline of this paper 

This is the fourth paper of this series and puts a new emphasis on regions in which the ambient pressure decreases with height 
like z"*. It also gives detailed analytical solutions of the dynamical problem for the first time. Papers II and III concentrated 
on why there are coUimated jets at all. Here we concentrate on the dynamic magnetic configurations generated. 

Paper I (Lynden-Bell & Boily 1994) showed that when an inner disk was rotated by 208 degrees relative to an outer disk 
field lines that had connected them splayed out to infinity in the absence of any confining external pressure. At greater angles 
there was no torque as the inner and outer disks were magnetically disconnected. 

Paper II (Lynden-Bell 1996) demonstrated that inclusion of a weak uniform external pressure led to a strong coUimation 
after many turns with the magnetic field creating towers with jet-like cores whose height grew with each turn. Again the inner 
disk was rotated rigidly relative to the outer disk. 

Paper III (Lynden-Bell 2003) was a refined version of a conference paper (Lynden-Bell 2001). In these the differential 
rotation of the accretion disk and external pressure variation with height were included and the shape of the magnetic cavity 
was calculated as a function of time. However the fields were not calculated in detail and the treatment used a static external 
pressure which was assumed to fall less rapidly than z^**. 

In this Paper IV we show that these quasi-static models are actually dynamically correct provided that the motions 
generated in the field lines by the twisting of the accretion disk never become relativistic, but we then explore the consequences 
of the magnetic cavity expanding into a region where the pressure variation approaches z^'^. This results in a dramatic 
acceleration of the top of the magnetic cavity along a cone whose angle gradually narrows at greater distances. The sudden 
expansion out to infinity when the twist exceeded a critical angle, found in paper I when there was no confining medium, is 
still present in modified form when beyond some height the external ambient pressure falls as z^'* or faster; once the field 
has penetrated to that region there is no longer a quasi-static configuration for the system to go to, so the jet accelerates to 
reach either dynamic ram-pressure balance whenever the background density falls less fast than z~^ or failing that relativistic 
speeds. 

1.4 Dynamics from statics 

In the standard MHD approximation the displacement current is neglected so curlB = 47rj. If the magnetic field dominates 
over any other pressure or inertial forces, then, neglecting those, the magnetic force density has nothing to oppose it, so 
j X B = and we deduce that the currents fiow along the lines of force. This implies that j =aB and as both divB and divj 
are zero a is constant along the lines of force. 
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We shall be considering problems with the normal component of magnetic field specified on an accretion disk at 2 = 
and the field does not penetrate the other boundary where an external pressure p{z) balances B^/Stt. The past motion of the 
accretion disk has produced a twist in the field-lines which emerge from the disk and return thereto further out. Those with past 
experience will know that the above conditions supplemented by expressions for B„ on the disk, the twists $ on each line and 
p[z), serve to define the problem, so the magnetic field is then determined. The other Maxwell equations curlE = —d'B/dct, 
and that giving divE, are not needed in the determination of the magnetic field. Notice that all the equations used to find the 
magnetic field do not involve the time. Thus if we specify the twist angles together with the normal field component on the 
disk and the boundary pressure p[z) at any time, the whole field configuration at that time is determined. Now let us suppose 
that we know how to solve this static problem but that Bn{R), ^{R) and p{z) are continuously specified as functions of time. 
Then the solution for the magnetic field in this dynamic problem is found by merely taking the sequence of static problems 
parameterised by t. Such a procedure will give us B(r, i) but in the dynamic problem the motions of the lines of magnetic 
force generate electric fields that must be found also. However these can be found almost as an afterthought because we can 
determine how the lines of force move and their motion is the cE x B/B"^ drift. Assuming perfect conductivity there is no 
component of E along B so E = B x u/c. Our knowledge of u on the accretion disk tells us how the lines move everywhere, 
which in turn tells us the electric field. 

We conclude that the crux of the problem lies not in difficult and dangerous dynamics but in the staid simplicity and 
safety of statics. That said we need whole sequences of static solutions that allow us to turn up the twists and parameterise 
the external pressures. Even the static problem is no walkover and we would have found it impossible to get general solutions 
were it not for the energy principle that the magnetic field adopts the configuration of minimum energy subject to the flux, 
twist and pressure constraints imposed at the boundary. 

Here we have already demonstrated why an evolution of the magnetic field structure through quasi-static models gives 
the solution to the dynamic problem. Section 2 details the specification of the relevant static problem and the methods of 
solving it. In Section 3 we solve it developing further the approximate method of paper HI. This gives us the mean fields 
within magnetic cavities whose shapes we calculate for any prescribed external pressure distribution. Emphasis is placed on 
solutions that access regions with p falling like z~* and the very fast expansions then generated. We then generalise our results 
to allow for a dynamic ram-pressure and discover the shapes of inertially confined jets. In Section 4 we calculate the detailed 
magnetic fields within the cavity. Section 5 finds the electric fields generated as the magnetic cavity grows and categorises the 
types of solution. Section 6 gives exact solutions of special cases with the dynamical electric field also calculated. 



2 THE MAGNETIC PROBLEMS TO BE SOLVED 

A magnetic flux P{Ri,0) rises out of an accretion disk on 2 = at radii up to Ri . The lines of magnetic force on a tube 
encircling the flux P eventually return to the disk at an outer radius Ro{P) after a total twist around the axis of <I>(-P). The 
magnetic field above the disc is force free with current flowing along the field lines and there is negligible gas pressure within 
the magnetic-field-dominated cavity. However the magnetic cavity is bounded by a surface at which an external pressure p{z) 
is specified. Later we shall consider the case of a dynamical pressure p{z,t). Our problem is to find the magnetic field and 
the shape of the cavity containing it when the functions P{R,0),${P) and p{z) are specified. Axial symmetry is assumed. 
We think of $ as due to the disk's past differential rotation and sometimes write <I>(P) — [Qd{Ri) — ^d{Ro)]t — Q{P)t where 
the suffix d refers to the rotation of the disk itself. The equation divB = implies that the magnetic field components in 
cylindrical polar coordinates {R, (j), 2) may be written in terms of the flux function P{R, 2) which gives the flux through a ring 
of radius R at height 2, and the gradient of the azimuthal coordinate 0, in the form, 

B = (27r)"VPx V0-f 5^0. (1) 

The force-free condition 47rj x B =curlB x B = then tells us that B^ takes the form, 

B4. = {2^R)-'P{P). (2) 

The function 13 is constant along each field line (so it is a function of P) and must be determined from the solution so that 
the total twist on that field line is ^(-P). Finally the azimuthal component of the force free condition yields the equation 

f) T r)P r)^ P 

V'P -V\nR\VP = Rf^i^^) + ^ = -P'iP)PiP) = -S^'RU- (3) 

This is to be solved within an unknown surface S given hy R = Rm{z) in which the field lies and on which B^ = 8np{z). 
The difficulties of this problem are: 

- the equation is non-linear 

- the function l3 in it is not known and can only be determined from ^{P) via the solution. 

- the bounding surface S is unknown. 
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Luckily there is a different way of tackling these magnetostatic problems. The energy of the magnetic field and of the 
external gas-pressure must be a minimum subject to the flux and twist conditions on the accretion disk. The pressure energy 
stored in making a cavity whose area at height z is A{z) against an external pressure p{z) is 



Wp= p{z)A{z)dz. (4) 

The energy principle that applies even outside axial symmetry is that W must be a minimum, where 



8ttW-- 



'B% + bI + Bl)Rd(t>dR + S>TTpA 



dz, (5) 



and B satisfies the flux and twist conditions on z = 0. Two exact theorems follow; they were proved in paper III by expanding 
a horizontal slice through the configuration first vertically and then horizontally. Earlier more specialised versions appeared 
in papers I and II. These theorems are true even without axial symmetry. Defining averages < .. > over a horizontal plane at 
height z, 



<Bl> + <Bl >=< Bl > +87r L{z) + [A{z)]-^ f A{z') 



)[dp/dz']dz' ] ■ (6) 



the final integral from z up is negative whenever pressure falls at height. Minimum energy for horizontal displacements of the 
slice gives, 

< B^ >= 8np{z) - {4n/A)dWo{z)/dz, (7) 

where 

4tvWo = / / BRB,R^d(j}dR, (8) 



and the integration is over the plane at height z. 

We shall show presently that after much twisting the magnetic configuration becomes very tall as compared to its width. 
As no extra radial fiux is created the radial field lines become widely spread out over the height of the resultant tower and 
the gradients with height likewise become small. When we neglect terms involving Br equation (|7J becomes 

<Bl>=STxp{z), (9) 

and with a similar neglect in the use of equation Q we find 

<Bl> /(2 - s) =< Bl >= 87rp(2), (10) 

where we have defined a dimensionless s (positive when pressure falls at greater height) by, 

s{z) = - / A{z){dp/dz)dz /[A{z)p{z)\. (11) 

J z 

Notice that s = for the constant pressure case and that (1 — s)/s is the constant d\n{A) / d\n{p) when A the cross sectional 
area varies as a power of p. 

Two possible approaches to solving this problem are: 

I. FORWARDS METHOD 

We use the energy variational principle choosing such simple trial functions that the boundary conditions can be applied 
and the resulting variational equations can be solved. We can then solve the problem for all choices of the prescribed functions 
P{R, 0), 'l'(P) and p{z) but the accuracy of the solution is limited by the imposed form of the trial function. 

II. BACKWARDS METHOD 

Here we solve the exact equation Q but make special choices of /3(P) and of the surface S so that we can solve the 
problem. Once we have the solution we discover the functions P{R, 0), ^(-P), and p[z) for which we have the solution. While 
this method is limited to a very few special cases, at least for them the solutions are exact. This enables us to check the 
accuracy and the validity of the more general solutions given by the Forwards method. 



3 SOLUTION BY VARIATIONAL PRINCIPLE 

This method was developed in paper HI and with an extremely crude but instructive trial function it was used in paper II. 
There we showed that if each field line turned N times around the axis and the magnetic cavity was taken as a cylinder 
of height Z and radius R then, if the total uprising poloidal fiux was F, very crude estimates of the field components 
are: Bz = 2F/{-kR?)]Br = F/ {\/2tv RZ) ; B^ = NF/{RZ). Squaring these estimates, adding the external pressure p and 
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multiplying by the volume -kR^Z we find Stt^W = F'^[{l/2 + N^-k'^)Z~'^ + [AR''^ + 8-K^pR'^F''^)Z] where p is assumed 
independent of z. Minimising over R gives R = {2tt'^p/F'^)~^''^ so R is determined by the external pressure and the flux. With 
this result the final round bracket in W reduces to 8-R~^ and minimising W over all Z we find AZ/R = \/l + 2N^tv'^ -^ ^/2-kN 
which shows a remarkable coUimation which improves with every turn! 

The method of paper III involved a much improved trial function which allows each field line, P, to attain whatever 
maximum height, Z{P) it likes, allows for the different total twists $(P) of the different field lines and properly accounts 
for the variation of external pressure with height. Because of the great heights of the magnetic towers generated after many 
turns, most of the field energy is high above the accretion disk and the detailed distribution of the flux in P{R, 0) no longer 
plays a part in the distant field. At each height, z, mean fields are defined by 

_ r^™ 

B^{z) = R-^ B^{R,z)dR, (12) 

Jo 

where Rm{z) is the radius of the magnetic cavity at height z, and A{z) = 7ri?^„. Also 

]b7\^A'^I \B,\2-KRdR = A-^j\dP/dR\dR = 2Pm/A. (13) 

The last expression arises because P is zero at both R = and R = i?^ and achieves its maximum Pm{z) at an intermediate 
point. We also define 

f =< Bf > /m\ (14) 

and 

J' =< Bl > Ib\- (15) 

although I and J are in principle functions of height, for tall towers we expect them to settle down to some typical values 
which we determine later. In what follows we neglect any variation of J with height; variation of I does not change the result. 
We now explain the basis of the trial function used. 

If a poloidal fiux dP is twisted once around the axis, it generates a toroidal fiux dP. If its twist is ^(P) the toroidal 
flux generated is (27r)~^<l>(P)dP. First consider distributing this toroidal flux uniformly over the height Z(P') to which this 
field line reaches as was done in paper III. Unlike the use of Z in the crude calculation in its new meaning it depends on 
the P of the field line. In a small height interval dz the element of fiux dP contributes a toroidal fiux \2-kZ{P)\~ ^{P)dPdz 
whenever the height is less than Z{P). The total toroidal flux through dz is contributed by all lines of force that reach above 
that height; that is by those with P < Pm{z), where Pm,{z) is the maximum value that P achieves at height z. Hence there 
is a toroidal fiux through the area Rmdz of 

R,nB4,dz = {2n)~^ / {$/Z)dPdz- 
Jo 
$(P) is specified in terms of the accretion disk's twist. Z{P) is to be varied. However while this makes a possible trial function 
it is not generally true that the fiux is so distributed. Indeed in the limiting case of the exact self-similar solutions discussed 
in the next paper the twist is strongly concentrated toward the top of each field line. The extreme alternative is found by 
putting all the twist at the top of each field line. Then the toroidal flux in any height increment dz depends on the twist 
of those field lines whose tops lie in dz so RmB^dz — {2Tr)^^${P,n){—dPm/dz)dz. However the truth must lie between the 
extremes of uniform twist and all twist at the top of each line. We therefore take the geometric mean of these two expressions 
for the toroidal flux distribution and obtain: 

fPm 

RmB4. = {2T^r\HPrr.){-dPm/dz) / {<^ / Z)dP]^'^ . (16) 



Now Z[P) is in the variational principle and Pm{z) has the property that Pm[Z{P)] = P, so P,„ is the functional inverse of 
Z{P); so whenever Z{P) is varied Pm{z) automatically varies in concert. We now omit the B^ term in the energy principle as 
it is much smaller than the others once the towers grow tall, cf the crude estimate above equation 11211 . Writing W in terms 
of the expressions discussed above 



(47r)~V^<E>(P,„)(-dP,„/dz) / ^/Z.dP + ApP'^/A + d,TTp{z)A 



dz (17) 



Varying A{z) gives 



A = ttP^ = IP^/yf2'Kp{z) (18) 

which using 11311 and 1141 1 is equivalent to < B^ >— Sirp, as found in equation ||5J. Putting this expression back into 
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equation(17) we find that the last two terms there combine to make PmdTl/dz where 11 — P 47-y/27rp(V)dz'. We now 
integrate both terms in W by parts and then rename the dummy variables in the integrals that result from the first term; 
these operations give 



F rF rF 



SttW- = (47r)- V' / [^{Pm)/Z{Pm)] j ^{P)dPdPm+j HdP^, (19) 

Jo J Pm Jo 

where IT is merely Il{z) re-expressed as a function of P,n 



n(P„) = / 4Iy^8npiz)dz (20) 



Most remarkably our adopted geometric mean between the uniform twist and top- twist cases has resulted in a H^ of precisely 
the form obtained in the uniform twist case of paper III except that the first term is now exactly half of its former value. Since 
this clearly reduces W the new trial function is clearly better than the old in that it brings us closer to the true minimum. 
Minimising W over all choices of Z{Pm) we find from 11911 and 12UII 



Z^^87Tp{Z) = [jV(167r/)]P^$(P^)' (21) 

where 

$' = <i>{Pm)P-' / $(P)dP (22) 

JPm 

Equation H21|l gives the function Z{Pm) and its inverse Pm{z) in terms of the given p{z) and <1?(-P). In line with section two's 
introduction we write ^ = fit and $ = fit. From equation (21) we notice that as a result of this t dependence $ and Z grow 
with t at each value of Pm- This growth is at constant velocity if the external pressure is independent of height. However if 
that pressure falls off with height the velocity dZ/dt accelerates, e.g. if p oc (a + Z)~^ with n = 2 then Z has nearly constant 
velocity until it reaches 'a' but at greater heights it behaves as t^ with constant acceleration. The shape of the magnetic 
cavity at each moment is given by plotting Z against Rm{Z). We find the relationship by substituting Pm from equation (18) 
into equation (21) to give the equation for Z as a function of Rm- Viz Z = [J/{'iV2I)]Rmfl[TvI~^Rm\/'2iTp{Z)].t, where fl is 
evaluated at the value of Pm indicated in its square bracket. Dividing equation (21) by equation (18) we find the coUimation 
at any given Pm is just 

Z/Rm = [J/(4V2/)].r!(P„)i, (23) 

which still grows linearly with time even when the external pressure decreases with height. However p{z) should not decrease 
too fast else equation 12H will not have a sensible solution. To understand this we see from equation 11221 that Pm$ decreases 
as Pm increases since the twist "I>(Pm) is greatest for the field lines rising nearest the centre of the disk. Also Pm decreases 

— 2 

at greater heights since less magnetic flux reaches there. Hence Pm^ increases with height Z. This must be true of the left 
hand side of equation I2H too and indeed it is obviously so when p is constant so there is then a sensible solution. However 
should p{Z) decrease faster than Z~^ the left hand side of l|21|l would decrease with Z so there would be no such solution. 
If the cavity accesses regions in which p decreases a little less fast than z~ then the field rapidly expands outwards and 
this gives a most interesting jet model. For example if p oc (a + z)~^ ~ ' then Z would grow with a very high power of t, 
t^' at each Pm so high expansion speeds would be achieved quite soon; however ram pressure will increase the effective 
pressure at the jet's head. This is explored in a later section. When relativistic speeds are achieved our approximation fails 
but analytical progress with relativistic jets can now be achieved following the lead of Prendergast (2005). For static pressures 
with p = po/[l + (^/fl)]" we illustrate the motion of the jet-head in Figure 1. For Z < a the velocity is almost constant but 
as the jet encounters the decrease in pressure it accelerates. For n — 2 the acceleration becomes uniform, but for larger n it 
increases and for n = 4 infinite speeds would occur in finite time if our equations still held. For n > 4 the graph turns back 
so the solutions indicate an infinite speed and then turn back giving no solutions for later times. These results are of course 
modified when ram-pressure is included as described later. 

3.1 CAVITY SHAPES FOR SPECIFIC MODELS 

We have the solution for any specified distributions of p{z) and for any twist $(P) — Q,{P)t. However before we can 
draw any cavity shapes we must also specify fl{P) or equivalently fl{P)- These two are connected through the definition 
P[0(P)]2 = n{P) Jp n{P')dP' = -{l/2)d/dP[Jp n{P')dP'f. This relationship is easily inverted. Evidently /J Q.{P')dP' = 
[/J 2PWdP]^''^ so finally fl{P) = Pff[Jp 2PTfdP]~^/^. We expect the behaviour of n(P) near P = P should be propor- 
tional to (F — P)^ since on the disk P^ = there. From the definition above this implies that Q, should be proportional to 
(P — P)'' near there. It is also true that if Q, tends to fio as P tends to zero then Q will be proportional to P~^' ^ near there. 
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Figure 1. Figure 1 shows tlie time evolution of the height of the jet in different pressure environments. On the left for the lowest graph 
the jot penetrates a constant pressure environment, n=0, at constant speed; the higher graphs are for n=2,3,4,6, and n=2 has almost 
constant velocity up to Z=a but feels the pressure decrease and thereafter proceeds at almost constant acceleration. At higher n the 
acceleration increases and there are no continuing solutions with magnetic pressure in balance beyond n=4. In fact the n=6 curve turns 
back unphysically to earlier times because the external pressure is too weak to withstand the magnetic field at larger heights. In reality 
the magnetic field springs outward at such a speed that dynamic ram-pressure comes into play. On the right we see the equivalent figure 
for ram-pressure with density profiles m=0,2,3,4,6,8 see equation 1261 



n 



-20 -10 10 20 
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Figure 2. The two figures on the left give the temporal evolution of the cavity shapes for the Dipole (left) and Simple models for the 
n=3 pressure distribution. The times corresponding to a given jet height are given in figure 1. Next we have the Simple model at a later 
time (notice the scale change). These are contrasted with the evolution of a magnetic cavity in a constant pressure environment n=0, on 
the right; here at each time the cross sectional area is always smaller at greater height. 
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Figure 3. < -B? > / < B^ > plotted as a function of \/z/Zfi for pressure-confined Simple-model jets with from the top n=0.2,0.5,l,2,3. 
Both field components become zero at the jet's head 

The simplest expression with these properties is the Simple Model Ti = (2)~^''^no(P/-F)~^'^^(l - P/F)^^^. This steadily 
falls from its central value and corresponds to f2(-P) = f2o(l — P/F). 

A model with a more physical motivation is a uniformly magnetised rapidly rotating star or black hole giving an initial 
flux distribution on the disk P{R, 0) = FP? /a^ for R < a and P = Fa/ R for R > a. This Dipole Model has a sudden 

reversal of field at i? = a. We combine this with a rotation Qd = ^* R < a;Qd ~ 0,t{R/a)^"^''^ R > a. Setting 

Y = P/F, fl = n4l- y^/^] and on integration n{P) = n^Y-^/'^y/{l - Y^/^)[3/5 -Y + (2/5)Y5/2]. 

The shapes of some of the dynamic magnetic cavities generated by combining these Q{P) distributions with simple 
ambient external pressure distributions p{z) are illustrated in Figure 2. 

Figure 2 shows the time evolution of the magnetic cavity for the n = 3 pressure distribution for the Dipole and the Simple 
models. They do not differ much. Very different cavities and velocities arise when the pressure distribution in the ambient 
medium is changed. In particular larger n albeit less than 4 gives a fatter jet at a given length and a longer jet at a given time 
as illustrated in figure 2, however the coUimation as determined by the length to width ratio is governed by equation 1231 
and so remains about the same at a given time. With the shapes of the cavities known it is now possible to calculate their 
areas at each height for each external pressure and thus discover how s{z) varies with height. This is interesting as, from IIUH . 
< B^ > / < -Bf >— 2 — s. For z » a the distribution of 2 — s with {z/Zhf is given for various values of n using the Simple 
Model. These graphs illustrate that the twist is less for the higher values of n but for them it is more concentrated toward 
the top of the jets. The integrations for equation lllH were performed by expressing the area as a function of pressure via 
equations (O, (HTJ and ll^. The final function for z » a which is plotted is 2 - s = (l/2)[(4-n)/(n-l)]a::''[(l-x'=)/(l-x'')] 
where x^ = {z/Zh);b = (4 - n)/3; c = 4(n - l)/3 



3.2 FAST JETS - RAM PRESSURE 

Our most suggestive finding thus far is that there are NO quasi-static solutions of large total twist $ when the external 



pressure falls like 



or faster. This result is easily understood; A purely radial field in a bottom-truncated cone r > a 



of solid angle u) will fall like r~^ . If the total poloidal flux both outwards and inwards is F the magnetic field would be 
2F/[ujr'^\ and it would deliver a pressure F'^ /[2-KUj'^r'^\ on the walls. The total field energy would be F^ /\/l-KUja] and an equal 
amount of work would be needed to inflate the magnetic cavity against the external pressure so the total energy would be 
F^ /[-KLija]. If instead we had a pure potential field with the same flux its field would fall as r~^ '^'^' so its energy would be 
about _F^/[27r(2Z + l)uja\ where I is the order of the Legendre polynomial that fits into the solid angle; so the energy of the 
purely radial field to infinity and back is only about 21 + 1 times the energy of the potential field. In practice I = {A-k/luY'^ . In 
paper V we find that a total twist of the upgoing flux relative to the downcoming flux, of n / sin[9m/2\ is sufficient for the flux 
to extend to infinity within a cone of semi-angle 6m- When the pressure falls with a power a little less negative than minus 
four we already demonstrated that continued twisting leads to the top of the tower or jet head advancing with a very high 
power of the time. However that was on the basis of a static pressure which will be enhanced by the dynamic ram-pressure 
once speeds comparable with the sound speed are achieved. Applying Bernoulli's equation in the frame with the jet-head at 
rest we have v^ /2 + [7/(7 — l)]p/p ^constant, so at the stagnation point ps — p[l + (7 — l)/2^.pv'^ /py'^'^^^' . If p >> pv^ we 
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get ps = p + pw^/2. Although that is much used we are more interested in the general case. When the flow becomes supersonic 
a stand-off shock develops so the above formula needs to be modified. Just behind the bow shock the pressure and density 
are given in terms of the upstream Mach number M by P2 = p(7 + l)~^[27Af^ ^ (7 ^ 1)1; p2 = p(7 + l)M^/[(7 — \)M^ + 2]; 
also P2V2 = pv. Applying Bernoulli's theorem after the shock results in a formula which for Mach numbers of two or more 
is well approximated by ps = r{'y)pv^. Here T = [(7 + i)/2](^+i)/(^-i)7~t/(^-i) which is 1,0.93,0.88 for 7 = 1,4/3,5/3, 
respectively. The full formula for ps has a further factor on the right [1 — (7 — 1)/(27M )]~^' 't~ ' which clearly tends to one 

at large M^. For 7 = 4/3 it ranges between 1.49 > 1.05 as M ranges from one to three. For 7 = 5/3 the range is from 

1.41 > 1.03. A useful global formula that somewhat approximates the trans-sonic pressure but is good in both the static 

and strongly supersonic limits is 

Ps=p + pv\ (24) 

This is the stagnation pressure that will be felt at the head of the jet. Away from the head the pressure is reduced both 
because the magnetic cavity expands less rapidly in the direction of its normal there and because a considerable fraction of 
the velocity is now parallel to the surface of the cavity. We return to this in the paragraph on Inertially Confined Jets. 

3.3 MOTION OF THE JET-HEAD 

From equations 12H & 1221 the position of the jet head at time T is Zh where 



Zl^ySTTps = J^[167rJ]"'n(0) / n{P)dP.t^ (25) 

Jo 

where, because we are dealing with the head, we have replaced the pressure with the stagnation pressure and set Pm = 0. With 

the stagnation pressure given by (I24II but Zh replacing v, equation 11251 is now an equation of motion for the jet-head. As only 

Zh occurs in this subsection we shall usually drop the suffix h which will be understood. When formulae are to be used in other 

sections we shall resurrect the suffix so that they can be lifted unchanged. Squaring 1211 p{Z) + p{Z)Z^ — L'^{t/Z)'^ where 

p{z) and p{z) are the undisturbed pressure and density at height z and L is the constant (l/2)(87r)~^" J^/~^r2(0) J Q{P)dP. 

We are interested in this equation when pressure falls with height. If the initial velocity given by neglecting the Z^ term is 

subsonic then initially we shall have results very similar to the quasi-static case illustrated on the left of figure 1 but all those 

solutions accelerate as the pressure decreases so the velocity will become sonic and then supersonic so that the p{Z)Z^ term 

will dominate over the falling pressure. Now neglecting the pressure and taking the square root we have Z'^[p{Z)Y'^Z = Lt^. 

Setting 

p{Z)^po[l + {Z/afr'"/', (26) 

so that the density behaves as 2"™ at large heights, we integrate to find 

[1 + (Z/a)^]'-'"/'' - 1 = (1 - m/6)Lii^ (27) 

where Li — Lp^ a~^. For all m the solutions for small heights are Z — a(Li)^''^f and for large heights 

Zh = a[(l - m/6)Li]''''^~'"'i'''''®-™' m < 6. (28) 

Once again such solutions accelerate but less rapidly; however if the density ever falls as fast as height to the minus six even 
the ram-pressure is too weak and the solutions rise formally to infinite speed before failing altogether when t — [6/(771 — 
6)Y''^F^'^L^ . The trouble arises because the rapid fall in density gives too small a ram-pressure to resist the acceleration 
of the jet. In practice either the excess inertia due to relativistic motion or the fact that the density does not fall below 
intergalactic values will avoid this behaviour. In the former case we need to develop the relativistic MHD jet theory. In the 
latter we may use our 771 < 6 model modified by replacing the density by the intergalactic one whenever our formula yields 
a smaller value, i.e. at Zh > a[(po/pO'^ ~ 1]^". As it reaches this region the velocity reaches Zh = [L{pi)~^''^]^'^ and 
thereafter it remains at that value. So Zh{t) is given by formula H28|l displayed in Figure 1 (right) until it reaches that region 
but then maintains its constant speed. 

3.4 INERTIALLY CONFINED SUPERSONIC JETS 

Landau and Lifshitz in their Fluid Mechanics book paragraph 115 give an elegant theory of supersonic flow past a pointed 
body, and the pressure on the body may be found by using the stress tensor p{d(j>/dR)^ of the velocity potential given in their 
formula 115.3. However there are two drawbacks. Firstly it is not likely that our jet will constitute a pointed body rather 
than a blunt one, even if it could be treated as a body at all, and secondly the resulting formulae involve integrals over the 
shape of the body which we can only discover AFTER the pressure is known. We shall circumvent such difficulties while 
maintaining momentum balance by treating the medium into which the jet penetrates as ionised dust each particle of which 
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collides with the magnetic cavity. Taking axes that move with the top of the field lines labelled by P i.e. with velocity Z[P) 
we find a ram- pressure on the cavity wall at height Z of 2pZ^cos^6 where cos^O — . ^. ^. . The factor two arises from the 
assumption of specular reflection from the cavity wall and at normal incidence the formula gives a factor two more than the 
stagnation pressure at the jet head. This suggests that a dead-cat-bounce off the cavity wall may be a better approximation 
than a specular reflection so we shall omit the factor two in what follows. However even the resulting formula is hard to apply 
in practice so we now proceed to simplify it further. Most of our jets are not far from parabolic at the front, in which case 
Zh — Z = K,{t)Rl^. Then {dZ/dRmY = 4(Zh — ZY / R^. Hence we shall adopt for the dynamic pressure at Z{P) 

Pi = pZ{Pf/[l + 4(Zh - Zf/Rl,] (29) 

To find the shape of a dynamically confined jet we need to solve equation I21I I with the dynamic pressure pd replacing p and 
use equations 12311 and 1261 1 for Rm and p{z). For Z » a equation 1211 1 becomes 



Z'-'"''2^^87rpoa™/[l -H4(Zh - Zy / Rl,] = [J^ / [I&t^ I)]PM^ t" . (30) 

Except near the jet-head Zh~ Z » Rm/'2 so we may neglect the one and , using (23) for Rm the above equation simplifies to 

Z^-"^/'^Z = KiQ^t^{Zh- Z). (31) 

For an initial orientation we take Zh >> Z. We may then integrate directly using our former result that Zh oc fi/('^~™/s) and 
obtain Z oc i(5-2m/3)/[(i-m/6)(4-m/2)i ^ ^'^^^^^ j^- foUows that Z grows with a higher power of t than Zh (at least for Z smaU) 
indeed X — Z/Zh oc fi/t''"'"'^). When Z is sizeable the above overestimates Z so it overestimates dlnX/dlnt which must in 
reality lie between zero and (4 — m/2)~^ which is itself less than (1/4) (1 — m/&)~^. We get a better estimate of the behaviour 
of X by writing equation I3H in the form 

^4-m/2j^ _ m/6)(l - X)Y^[1 + (1 - ■m/6)dlnX/dlnt] = K2t, (32) 

where at given P, K2 is constant. Now the second square bracket above only varies between one and 5/4. If at lowest order 
we neglect its variation, we may solve for t{X). We may then evaluate 

dlnX/dlnt = (1 - X)/[4(l - m/8)(l - X) + X], (33) 

Ij32|l then determines t as a function of X. Evidently the ti/t^^'"/^) behaviour of X persists approximately until X is near 
unity; however Zh — Z oc ^^'('^"'"'^'(l — X) and this actually grows because with X near unity the second square bracket in 
equation ^ is one and so Zh - Z (x [K2{1 - m/6)]-'^^^'^-'"^'^'^ X^^'""^^^'^^"'"'^^ {I -- X)-"'^^^-"'\ so as X grows Zh-Z, the 
distance from the head actually increases. Z never catches up with Zh despite the fact that X tends to one. All the above rests 
on the premise that Zh — Z » Rm/2 so we now investigate the behaviour when Z is close to Zh so that X is close to one 
but not so close that (1 — X)Z/Rm has to be small. Writing equation l)3()|l in terms of X, but omitting the dlnX/dlnt term as 
this vanishes with X near unity, and dividing it by the same equation for Zh we find 1 + 4:{Zh — Z) /Rm = X ~"^ /[K4{Pm)] 
where K4 is Pm^'^ (Pm) divided by its non-zero value when Pm tends to zero. At constant Pm, {Zh — Z)/Rm clearly increases 
as X increases but it tends to the limiting value (l/2)y'K^ — 1 as X tends to one. As K4 depends only on Pm the shape of 
the jet-head is determined by the n(P) function though its size and position depend on t also. Thus the whole head of the 
jet and all parts with X near one grow self similarly. The shape of the these parts is even independent of the details of the 
density fall-ofi^ embodied in m but the size of these parts of the jet is proportional to t"^/(^-"^) go the rate of self-similar growth 
depends on m. As explained above the whole length of the jet grows with a different power of the time so it is just the part 
with X close to one that grows self-similarly. To find the shape of the whole cavity we notice that elimination of Omega-bar 
between equations 12311 and 1301 1 or 12111 coupled with use of 13311 leads to an expression for Pm in terms of Z, Zh, Rm, t 



^"^ ^ 21(1- ^) '-'^'"^Z*^ 



^^ (l-m/6)(l-2-/Z,0 



[A{l-m/8)(l-Z/Zh) + Z/Zh] 



87rpoa'"(Z3+a3)-'"/3 

[i + A{Zh^z)yRi,] ^■^ > 



This value of Pm is substituted into the Omega-bar of equation (23) to yield the equation relating Z and Rm at each time, 
Zh{t) being already known. Thus we get the shapes of the inertially confined jet cavities. Figure 4 displays one of these at 
two different times. 



4 THE FIELD IN THE MAGNETIC CAVITY 

Our solution of the variational principle has improved on paper III and given us dynamical solutions for the cavity's shape 
and the mean field at each height. We now seek the detailed field structure within the cavity. At each height z we define A to 
be [R/Rm{z)] so P may then be written 

PiR,z)^Pm{z)f{\). (35) 
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Figure 4. An inertially confined jet cavity at two diff'erent times, the inner edge of the black gives the cavity in an ni=3 density 
distribution for the Simple model magnetic flux. The times to these heights are given by equation I27i 



Since by its definition Pm(z) is the maximum value that P takes at height z, it foUows that / achieves its maximum of unity at 
each height. Furthermore, as P is zero both on axis and at the surface of the magnetic cavity, it foUows that /(O) = = /(I)- 
Thus /, which is positive, is highly circumscribed rising from zero to one and falling again to zero at one. Although / may 
in principle depend on height (especially near the disk or at the top), nevertheless so circumscribed a function is unlikely to 
have a strong height dependence. We shall make a second approximation that an average profile will do well enough over at 
least a local region of the tower's height. Thus we adopt the form given in equation 13511 with / a function of lambda but 
not of height. Now in the tall tower limit both Rm{z) and Pm(z) are only weakly dependent on z, so squares of their first 
derivatives and their second derivatives may be neglected. Equation @ then takes the form 



4P„i?-'A.a7/9A" = -P.dP/dP = -(9/?79A)/(2P™/ ). 



(36) 



Multiplying by 2P™/' and integrating dX we find /3^ = 4P^R-^{J^ /'^dX - A/'^), so /?(P) is of the form /?[P™/(A)] = 
(2Pm / Rm)G{X) . Evidently /3 is a product of a function of z and a function of A but it is also a function of such a product. It 
is readily seen that a power law form for /? achieves this and we readily prove that the power law is the only possibility that 
allows it. Hence we may write 



(3 = CiP" 



(37) 



Thus P' (3 ex pP^'^ ^ . Inserting this into equation (1361 and calling the value of A where / achieves its maximum of one, Ai, we 
find on separating the variables A and z both 



Xd^f/dX^ 



-C uf 



(38) 



and 
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Rnr = CiP^-", (39) 

where 

C^ = C?C|/4. (40) 

Now Pm (z) decreases because not all flux reaches up to great heights, so equation 13911 tells us that the radius of the magnetic 
cavity decreases there when i/ < 1 and increases when i/ > 1. Oi course this only holds once the tower is tall so that its lateral 
confinement is due to the ambient pressure rather than the flux profile on the disk itself. The constant C is determined by 
the requirement that / = 1 at its maximum because /(A) already has to obey the boundary conditions at and 1. Sometimes 
we find it convenient to use a = 2;/ — 1 in place of u. We already know how to calculate the shape of the magnetic cavity. 
Knowing A{z) and p{z) we can calculate s — —{Ap)~^ J A{z'){dp/dz')dz' at every height z. We now show that associated 
with each value of s there is a profile /(A). Equation 13811 governs the possible profiles. Multiplying it by —//A and integrating 
by parts we find J^iffdX = C^v J^ f'"\~^d\. However P™./' = A{2TTR)-^dP/dR = AB^ so the first integral is related to 
< B^ > and the right hand one is similarly related to < B^ >. Multiplying both sides by Pm/A we find 

< Bl >=u<Bl> (41) 

so comparing this with equation IIUH wc find s = (2^ — l)/i^, so 

i/ = l/(2-s). (42) 

Thus for each height we have an s and hence a v for that height and the associated profile is given by solving equation 
13811 with / = at A = and A = 1 and the value of C is determined by the condition that / is one at its maximum. Notice 
that by this means we have determined the weak variation of the profile with height (see Figure 3) . 

4.1 PROFILE FOR CONSTANT EXTERNAL PRESSURE 

From equations 1181 and I39I I constant p corresponds to u = 1/2, q = 0, a case of great simplicity. Integrating equation 
(38) with the boundary conditions that / = at zero and one, we find /' — — C'^(l/2)ln(A/Ai) and / = — C^(l/2)AlnA with 
Ai — 1/e. The requirement that / = 1 at maximum then gives C^ = 2e so 

/ = -eXlnX. (43) 

With / known we may now evaluate the dimensionless integrals / andj, 
I^ =< Bl > /Jb7\' = (1/4) J^{df/dXfdX = eV4, hence 

7 = e/2 = 1.359. (44) 

Likewise 



J 



fhf/X)dX , 

•'o - ^ = 0.798. (45) 



\ /;(/^/VA)dA 



Thus once we specify 0.{P) and p, our solution for the shape of the magnetic cavity Rm{z) is known via equations 12111 . 12211 . 
and 1231 1. and within that cavity the structure of the field is given via the poloidal and toroidal flux functions P and fH. 

P = P„(^)eAln(l/A) (46) 

/?(P) = 4(27r^p)'''*P'''2 (47) 

The magnetic fields are given by Br = -{2TvRy^dP/dz, B4, = (27rP)"^/3, B^ = {2nRy^dP/dR. These are surprisingly 
interesting, 

which is positive for A < 1/e, zero at A = 1/e, and negative for A > 1/e reaching the value — ePm/(7rP^) at the boundary 
A = 1. By equation H18^ the magnetic pressure precisely balances the external pressure there. Unexpectedly we find that Bz 
is infinite on the axis 

B, oc ln(P„/P), (49) 

however the flux near the axis is small because P behaves like R ln{Rm/R) there. Likewise the contribution to the energy 
from the magnetic energy-density near the axis behaves as P^[ln(Pm/P)]^so the infinity in the magnetic field appears to be 
harmless. A greater surprise comes from the behaviour of B^ near the axis. 
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2eR, 



B^ = =^^HR^/R), (50) 

whereas we expected to find S^ to be zero on axis, it is actually infinite! However the ratio 

B^/B, ^ MRl/R^)]-'/^ -. (51) 

thus although the field lines near the axis do wind around it helically, the helix gets more and more elongated as the axis 
is approached so the axis itself is a line of force. Faced with this example it is evident that the normal boundary condition 
that on axis B^ should be zero should be replaced by the condition that B^/B^ should be zero on axis. Since the field is 
force-free the currents flow along the field lines and 47rj — /3'(P)B = (Ci/2)P~^'^B. As P is zero on axis j is even more 
singular on axis than the magnetic field. Nevertheless the total current parallel to the axis and crossing any small area is 
finite and tends to zero as the area shrinks onto the axis. The exact Dunce's Cap model of section 6 gives the same field 
structure close to the axis. In the next section we find that the infinite fields on axis are replaced by large finite ones when 
the pressure decreases at greater heights. Current sheets are of course a common feature of idealised MHD and we have one 
of necessity on the boundary of the magnetic cavity. The infinite field strengths and current densities on axis encountered in 
this solution suggest that very large current densities occur close to the axis in reality. The lack of sufficient charge carriers 
there will lead to a breakdown of the perfect conductivity approximation with large EMFs appearing up the axis resulting in 
particle acceleration along the axis. What observers 'see' as a jet may be just this very high current-density region where the 
particles are accelerated, i.e. only the central column of the whole magnetic cavity which may be considerably wider, cf the 
observations of Herculis A, Gizani & Leahy (2003). 

4.2 PROFILES WHEN PRESSURE DECREASES WITH HEIGHT 

Even when the pressure varies we know how to calculate the shape of the magnetic cavity, so all we need is the profile of the 
field across the cavity. So we need to solve equation 1381 . For v = 1/2 equation 14311 gives the solution, but we need it for 
more general i^. While this is easily computed an analytical approximation is more useful. For ;/ < 1 a good approximation 
is given by noting that for a — 2i/ — 1 small — InA ~ (1 — A")/q; using this in our a = solution for / suggests the form 
/ oc A(l — A")/o but a better approximation is given by / — g{X)/g{Xi) where 

5(A) = Q-^A(l-A")/(l + a2A"), (52) 

and Ai is given by g'{Xi) = 0. Notice that the denominator in g is constant when a = so it makes no difference to 
that solution. This form for / automatically satisfies both the boundary conditions and the one-at-maximum condition. The 
equation that gives Ai is a quadratic in A" 

1- [a + l-a2(l-Q)]A?-a2Af =0. (53) 

We use it to find a^ in terms of Ai which we determine below 

a2 = [1 - A?(l + a)]/[A?[l - (1 - a)A?]]. (54) 

The logarithmic infinities in the fields on axis found in equations 14811 - 150(1 occur because of the log term in equation (46). 
The form of equation 15211 shows that large finite fields replace those infinities when a > and so B^ -^- on axis. 
For a = 1 there is another exact solution in terms of the Bessel function Ji , 

f = kiVxJi{k^Vx)/[koJi{ko)]; (55) 

here feo ~ 2.405 is the first zero of the Bessel function Jo and fci — 3.832 is the first zero of Ji. The value of Ai is therefore 
[2.405/3.832]^ — 0.3939. This is not too far from the value l/e=0.3679 obtained for our a = solution. This suggests that 
linear interpolation i.e. Ai = e~^(l + 0.07073a) and a2 determined via equation 15311 will give a good fit and indeed numerical 
computations show the fit is excellent all the way from to l,the Bessel case. In the latter the expression fcoJi(fco) = 1.249. 

The integrals / and J can be evaluated for the Bessel solutions and we give their values for the a = case in brackets 
for comparison. / — 1.179(1.359); J = 1.098(0.799). As expected / and J do change with a but remain within 15% of their 
means. 

Solutions to equation 1381 for large a — 2v ~ 1 were used in paper I to solve a different problem but that method works 
well and can be extended to work for all a above unity. Approximate solutions for large v of the form / = 1 — v~^ ln(cosh A) 
where A = 2Ai(A-Ai)- Aln(cosh[i/C(A-Ai)]) with A = (6Ai)-i; Ai ^ch'^e" ~u; C ~2X\^^; Ai ~ f are deduced in the 
appendix as well as the generalised form 

/ = 1 - H'^ In(coshA); ....; A = C*(A - Ai) - Aln(coshA), (56) 

with H ^ u ~ 1/42 - 5/26. !/"^ h = VITu- A = 1/2. i/"^ - 1/9. i/"^ + 1/12. !/-^ Ai = 1/2 - 7/31. i/"^ + 4/41. i/'^ and 
C = hCX^ = 2iy + \n4 — 4/15m~^ . Notice that with A small equation 1561 with H — v Ct = v, CX^ reduces to 
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Figure 5. Computed profiles /(A) are plotted for i/ = 1, 3, 5, ... . A triangle is the limit v -^ oo. At the bottom are plotted the errors in 
our analytic approximate solutions derived in the appendix. Only for i^ = 1 do the errors rise above 2% and for v = 1 itself we have the 
exact Bessel solution. 

our former solution. Inclusion of the A term allows for the asymmetry around the maximum caused by the variation of the 
initial A in equation 13811 . These freedoms allow us to fit exactly, not only the curvature at the maximum and the boundary 
conditions, but also the gradients at which the solution reaches zero at the ends of the range. Some details are given in the 
Appendix. The solutions for / are plotted in Figure 5. By taking u to be given by equation 14211 in terms of s[z) we now have 
a solution of the form P = Pmf{X, v) with A = (R/Rm)^ and Pm, Rm and f all depending weakly on z, see equations IIBH . 
1211 1 and 1421 1. With P so determined the function f3{P) is found from $(P) by integrating along field lines. These are given 
by 4^ = 4^ = #. Hence 



dR 



Rd(j> 



dz 



(57) 



-dP/dz I3{P) dP/dR 

The equality of the first and last terms merely tells us that P is constant along field lines. The equality of the second and last 
terms gives us on integration with P held fixed 

2P$(P) 



P{P)- 



J{dlnf/dlnX)-^dz' 



(58) 



where the integration is along the curve P constant from one foot point to the other. While the above procedure is the one to 
use when P{R, 0), ^{P) and p{z) are specified there are special cases that are simpler. If we ask that / be strictly independent 
of z rather than that being an approximation, then equation l|89^ must hold exactly. Combining that with equation 1)18^ we 
would have p = 7^(27r)~^(7rC2A"°)"^/*^""' oc P^; from (11) this formula leads directly to s = 2a/{a + 1) and equation ^ 
gives us the simple relationship /3(P) — CiP^"^^'''^ in place of equation I58L nevertheless these equations with s constant 
imply that p vanishes, or becomes infinite, when the area A vanishes, which will not be true. The exactly separable case is 
too restrictive near the top when the pressure varies. 



5 ELECTRIC FIELDS, PUSHERS, FLOATERS AND SQUIRMERS 

We now determine the electric fields that occur when the accretion disk is in differential rotation. The velocity u of a magnetic 
field line is cE x B/B'^ and there is no E along B because of the perfect conductivity, hence E = — u x B/c . To use this 
formula we first calculate the velocity of our field lines. On the accretion disk this velocity is that of the disk itself so the field 
line which initially intersected the disk at azimuth i/) now has its outer intersection at Ro{P) at azimuth (l>d ~ "ip + Qd{Ro)t. 
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We now look for the angular velocity of the point at which this field line intersects a z =const plane. At each instant the field 
line obeys equation 15711 so on integrating the second equality there 

<f) = i) + nd{Ro)t + VL{P)tq{z.P) (59) 

where we have written Q{P)t for $(P) and q — [^{d[nf/dln\)~^dz/J{dlnf/d[nX)~^dz. Both integrations dz are to be 
performed with A varying so as to keep P constant. The final integral is to be performed from foot-point to foot-point. 
The physical meaning of q is the fraction of the total twist on the field line P that occurs by height z. Since the field line 
reaches a maximum height Z{P) and then returns, q will be double valued as a function of z. To avoid this it is better to 
convert those 2— integrations into integrations over A, in which q is single valued, rather than z. Such a conversion yields 
<1 = J^° Z'{P/f){f^)~^d\/ J ° Z'{P/f){fX)~^dX, where we have written P/f for Pm to demonstrate that the A dependence 
arises through / when P is held fixed. The angular rotation rate of the intersection of this field line with a z =const plane is 

<j> ^ Qa + n{P)q + n{P)tq (60) 

q arises from two causes i) any explicit dependence of Z' on t that does not cancel between numerator and denominator in q, ii) 
the change in the lower end point A of the numerator's integration. At constant P and z we have din Pm/dt = — (dln//dA)A, 
which gives us the rate of change of the lower limit of the numerator. When P{z) is either a power or constant the explicit 
dependence of Z' on i cancels between numerator and denominator so there is then no contribution from i). Of course had 
we left q as a z-integration there would have been no contribution from the end point but then the contribution from the 
t-dependence of the A in the integrand must be re-expressed as a function of P and z via /(A) — P/Pm{z,t), is unwieldy. 
Now R<j) at constant z is not the (^-component of the velocity of the line of force but the velocity of the point of intersection 
of the line with the z=const plane. When the line of force is only slightly inclined to the plane the difference is obvious since 
the velocity of the line is perpendicular to the line. A little thought shows that in the intersection velocity the component 
of the velocity of the line in the z-plane parallel to the projection of the field into the plane is exaggerated by the factor 
B^ / B^ relative to that component of the field line's velocity. Thus writing hats for unit vectors and b for the unit vector 
{BrR, + B^ij))/^/Bj^ + B? in the a-plane considered and u_l for the component of the field line's velocity in that plane, the 
velocity of the intersection is ux + [{B^ /B1) — l](ux.b)b. The component of this intersection's velocity along (j) is 

Ri = u4,{l + Bl/Bl)+UR{BRB4,/Bl). (61) 

Now by Faraday's law, the rate of change of the magnetic fiux through a circle about the axis in the plane considered gives 
the EMF so 

2-kRE^ = -dP/dct = -P/c = ~{u,Br - urB^)/c (62) 

This is a second equation for u and the third is just u.B = 0. Eliminating first Uz and then ur, we find u^f, — [{Br+ B'^) / B'^]R4>+ 
[BRBrf,/{B^ Bz)]P. The other components are readily found from equations 16H & 16211 . The remaining components of E follow 
from cE = B X u. Like the electric fields on the accretion disk these fields are not far from cylindrically radial and directed 
away from the surface of greatest fiux, P, at each height. 

We now turn our attention to categorising the different types of solutions. A Helium balloon needs a tether in tension 
if it is not to rise further. A water tank needs a support in compression if it is not to fall. If we cut a magnetic tower by 
a horizontal plane there is net tension across that plane due to the magnetic stresses if [{B^ — Br — B'^)dA/{8iT) > and 
net compression if that quantity is negative. We call these floaters or pushers respectively if they satisfy those criteria low 
down the tower. Floaters are are held down by magnetic tension; pushers are supported from the bottom by a net magnetic 
pressure. In the tall tower approximation the above criterion simplifies to (< 5^ > — < B^ >) > for a floater. Comparing 
this criterion with equations HUH and 1421 we see floaters have s>l,i^>l,a>l and cavities whose radii increase at greater 
heights while pushers have < s ^ 1, 0.5 < v !^ 1 and < a ^ 1. Their cavities' radii decrease with height. In a constant 
pressure medium all the solutions are pushers. If the pressure in a long narrow column supporting a weight is greater than 
the ambient pressure in the medium surrounding the column then any lateral bend bowing the column will be exaggerated 
as in the Euler strut problem. All our magnetic towers are stable against that bowing as their net support is atmospheric or 
less {< B% > + < BI > ~ < B'^ >)/(87r) -p = A'^ J A{dp/dz)dz < 0. This criterion corresponds to a > 0, s ^ 0. Thus 
squirmers do not occur in the magnetostatics of our systems which have pressure that decreases with height. 



6 EXACT SOLUTIONS, BACKWARDS METHOD 

Here we postulate the forms of l3{P) and of the bounding surface S for which we can solve the differential equation Q. 
After we have solved it we discover what problem we have solved by finding p{z) on the bounding surface and the twist 
angles ^{P) and the boundary fiux P{R, 0). The differential operator on the left of equation ^ separates in oblate or prolate 
spheroidal and rotational parabolic coordinates as well as in cylindrical and spherical polars. It is therefore important to take 
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the surface S to be one of those coordinate surfaces. Even after choosing one of those the problem is still non-linear; but there 
are interesting exceptions. We have already seen that the case /3'/3 =const. corresponds to a constant external pressure and 
it gives a linear but inhomogeneous equation. A second simple case is P{P) oc P oc /3'/3. This linear case is well known and 
much explored but is somewhat less interesting than the (3 ex P^'^ case. The combination /3'/3 oc P + Po though soluble has 
the idiosyncrasies of both the former cases without obvious advantages. A different approach with more interesting solutions 
is to accept the non-linear behaviour of equation Q and look for self-similar solutions. Most of the above methods can be 
generalised using perturbation theory. For example if the bounding surface is a 'cylinder' of slowly varying radius Rm.{z) one 
can solve the problem in terms of the scaled radius R/Rm as in section 4. 



6.1 SELF-SIMILAR SOLUTIONS 

Writing equation @ in spherical polar coordinates (r, 9, 0) and setting /i = cos 9 we find 
r^d'^P/dr^ + (1 - M^)a^P/c»M^ = ~r^P'f3{P). 



(63) 



Self-similar solutions only occur when P ex r~', I3{P) = CiP". 

If we now set P = r~^M{fj,) then the left hand side is r"' times a function of fi so /3'/9 which is a function of such a 
product, must equal such a product. Hence (3' (3 = (1 -I- l/VjC'lP^'^'^'^ . Thus with v = l-\- l/l the r~' cancels out and writing 
/ — M/M\ , where Af i is the maximum value of M 

l{l + 1)/ + (1 - )i')d'f/df,^ = -(1 + l/l)C',f+'^\ (64) 

This is a version of the equation of Lynden-Bell and Boily (paper 1) for whom the case with I small was of especial interest. 

To make contact with the work of section 4 we remark that for narrow jets it is only necessary to consider fi close to one. 
Supposing the walls at fi = Hm, we write A — (1 — /i)/(l — fJ^m), where both (1 — Hm) and (1 — /i) are small, the equation takes 
the form 



\d^f/d\' + l{l + 1)(1 - Mm)//2 = -~C\l + 1/0/' + '/'. 



(65) 



We have approximated 1 — At' as 2(1 — /i), which is good to 1 per cent or better for opening angles of less than 23 degrees. When 
the (1 — fim)l(l + l)//2 is neglected this equation will be recognised as equation 1)3 8|l of section 3. The particular solutions of 
equation 164^ to which we now turn are so simple that we obtain them without making this narrow cone approximation. 



6.2 THE DUNCE'S CAP MODEL 

Unexpectedly we shall use spherical polar coordinates centred at the top of the tower-like magnetic cavity at (0,0, Z). This 
use of Z is similar to the very crude description as a cylindrical tower in the introduction but now the tower will be conical 
pointing down at the accretion disk. This Z corresponds to Zh in sections 3 and 4. We measure 9 from the downward axis 
pointing toward the centre of the accretion disk (see Figure 6). We solve the force free equations within the downward opening 
cone cos6' ^ /irn which forms a dunce's cap configuration over the accretion disk. We take I — —2 which corresponds to 
our P (X P ' behaviour appropriate for a constant pressure. Equation l|64|l becomes 2f + {1 — fi )d f/dfi — —Ca,. Writing 
f = {1 — fi )m our equation becomes 

d/d/x[(l — H ) dm/dfi] — —C4, so, using the boundary condition that {1 — fi )m' is non-singular at /i — 1, we find, 
m — C4 J {I + ii)~^ {1 — fi)~^dfj,. On integration by parts we get 



P = r-{l - 11 )m = C4r-{1 - /i) 



3(1 + ^" 



^(^+'^)+il^l+M: 



1+ ^J. 1 — ^I^^ 



1 



M 



(66) 



Evidently P 
becomes P 



on fj. = Urn and on the axis. For tall narrow cones 6 is small, /i, jim are near one and the flux function 



\C4,r'^9'^[\n{9l^/9'^)]. At r = Z the maximum value of P is P so d = 'ieF/RJ where Ra = 9,nZ. Also 
/3(P) = ^2|C4iP'/' = C4r6l[ln(6'™/6»)]i/^ Also, see Figure 7, 

27rB, = {r^ey^dP/d9 = C4ln{9,^/9) - 1]; 2TvBg = Ci9ln{9,n/9); 2ivB^ = C4ln{9m/9)]^^^ (67) 

The constant d is determined by the external pressure C4 = 87r[7rp]^". Near the axis these fields have precisely the 
same behaviour as those we found in equations 1491 and 1501 . However the fields are in one sense even more interesting; 
they are all independent of r, the distance from the origin. By construction Bg,B^ vanish at the edge of the cone but Br is 
constant along its generators all the way up to the origin. The origin itself is then a point at which the infinite field coming 
up the axis turns back and splays out down these generators at constant field strength. This is the way the constant external 
pressure is opposed. With the origin such a remarkable point, some may suspect that the forces do not balance there, that 
the whole configuration might be held up by some sky-hook pulling at the apex. However P varies as r'^ and the Maxwell 
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Figure 6. The poloidal lines of force in the Dunce's cap model at two different times. The top advances at constant speed while the 
bottom only revolves in differential rotation 




Figure 7. Br, Bg and B^ are all independent of r inside the Dunce's Cap so we plot them as functions of 6/6max- While Bg is zero at 
both boundaries B^ is zero only at the outer one. Br passes through zero at an internal point. Sec text for the surprising behaviour near 

e = o 



stresses integrated over horizontal cuts through the cone give a net force that vanishes like r^ as the origin is approached. 
Thus no sky hook is necessary to balance forces at the origin. The radial fields reverse at an intermediate value 9i = 9m/\/e 
(for narrow cones) and some may suspect reconnection there but the B^ field component maximises at that value of theta so 
there is plenty of magnetic pressure to prevent reconnection. With Hm close to one the radii on the accretion disk are R = 6Z 
so the flux coming up within a circle of radius R is P{R,0) = Fe(i?^/i?^)[ln(i?^/i?^)]. The twist function ^{P) is given by 
integration along field lines which is simple for this self-similar field. rdO/Bg = rs'ai6d(j)/B^ yields on writing 77 = 6m/0 
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$(P) = n{P)t = [Z/Ri) / (ln7?)-i/"dr? ~ {Z/Rd)[{T)^ - r,o){\ni^^)-^'\ (68) 



where the integral is evaluated between foot points at Ri^o{P) of the field line and r]o^i{P) = Rd/Ro,i{P)- Therefore 770,1 are 
the roots for rj of 

P/F = 2eri'^ In r; (69) 

When the inner foot is much further in than the outer foot then Qi » Qo so then n(P) is close to Qi{P). This allows us to 
estimate the angular velocity of the accretion disk required. In the central parts it behaves as R~^[ln{Rd/R)]~^'^, so the circular 
velocity of the disk is almost constant but drops to zero at the very centre like [ln(_Rd/^)]~^ • From Equation 1691 both roots 
rjo^i are independent of time. Hence from equation 16811 Z ex f as expected from section 4 so we may write Z = Vt. Writing 
our solution for P in terms of cylindrical polar coordinates now centred on the centre of the disk, P = eF{R? /R^) \n{R^/R^), 
where Rm{z) = -Rd(l — z/Z) for z < Z and it is just the dependence of Z on time that generates the evolution of the system. 
(3{P) = VSeFP^'^ /Rd with ^{P) given by equation 16911 . The magnetic fields are given by equation 1671 1 so we now calculate 
the electric fields following the method of section 5. The simplicity of the expression for P means that q is simpler and in place 
of equation ^ we have (/) = V + f^d(-Ro)t + (Vt/Rd) J^'"'''^ (In r7)'^''^dr;. Now i?„ = (1 - i^)Rd so R„, = zRd/{Vt'^)4:. At 
fixed height a point only remains on the same field line if _R = ~P /{dP/dR) = —P /{2-kRBz) Using these and difi^erentiating 
S at constant z and P 



Rm/R 



z PVt 



^-nd + VRd^J^^ (ln,)-/^d,+ [ln(i?Wi^)]-^/^ [f,+^^^R-^^) ■ (70) 

As before Equations 1611 and 1621 give the velocities of the lines of force as, 

^0 = [{Bl + B'i)/B^]R^+[BRB4,/{B''B.)]P 

UR = [bI/{BrB^)][R4>^u^{1 + bI/bI)] (71) 

u, = {B,/Br)ur + P/Br 

Also cE^ = ~'P/{2nR) while the other components follow from cE = B x u. Notice that the meridional field velocity ua/ — 
{ur,Uz) is not perpendicular to the vector meridional component of the magnetic field. It is the full vector velocity that is 
perpendicular to the complete magnetic field. We have explored the Dunce's Cap model in all this detail both to demonstrate 
its interesting field structure and to show how we get the fields in an explicit example. Figure 5 gives the poloidal structure 
of the magnetic fields which wind around these surfaces of constant fiux. 



7 CONCLUSIONS 

It is very remarkable that these simple analytic methods allow the solution of these highly non-linear problems with variable 
domain, not just for a few special cases but for all twist functions Omega and all pressure distributions. Furthermore we 
get time-dependent solutions for all time. The main secret lies in the variational principle coupled with a good form of trial 
function. Once the time-dependence becomes relativistic we lose this tool so the problems will become harder. However when 
the relativistic motion is only important for the rotation it is likely that another variational principle may exist. An important 
problem is to seek it out. In the foregoing we have explored the simplest case in which the magnetic cavity is empty and 
the Poynting fiux carries both the energy and the momentum of the jet. I believe that the results justify the assertion that 
this simplest case has remarkable similarities to the observed world and this suggests that the extra complication of material 
winds is not essential to explaining the main phenomena. However especially for Pulsars Relativistic rotation is an essential 
ingredient of the problem not covered here. 
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APPENDIX A: A USEFUL MATHEMATICAL TECHNIQUE 

We wish to solve a highly non-linear differential equation such as 

(1 - m')0 = -1(1 + 1)/ - c^uf-' ; ^ = 1 + } (Ai) 

under the boundary conditions that / is zero at given end points and that C is so chosen that the maximum of / is one. 
There should be only only maximum of / between the end points. We rewrite the equation in the form 

where 

S^=l{l + l)(l-f)+C''{l-f''). (A3) 

We shall be particularly interested in the strongly non-linear case with u large {I small) so we start with u large but later 
extend the technique down to i/ « 1. As we shall come across several variants of essentially the same problem we describe the 
technique in terms of the more general problem in which IX2I is replaced by 

d I df \ d r„, j.-.,2 



with S^{f) zero when / reaches its maximum of 1. We shall assume that S{f) depends on parameters such at u or C in lAH 
and that we are primarily interested in the solution when u is large in which case dS^ /df is strongly peaked close to / = 1 
and away from there it is small. Under these conditions (which hold for IjAlfl ') most of the deceleration —d f/dfi of / occurs 
in the region in which / is close to one. Let that maximum he a.t /j, — fii. Then at /ii 

say (A5) 



(dV/dM')i = ^ 



df 
where Gi — G{fii). Near p = /xi 

(A6) 



(A7) 
(A8) 

(A9) 

and G'l = (dG/dn) at ni. 

Writing <A8l l as a cubic leads to simpler mathematics later. Near fj, = fj,i where most of the 'deceleration' of / takes place we 

see from 1A4II that 

d fdfY _ 2SdS/dfM ^^^^^ 



f = i- 


- 29{fJ-- tJ-if 










and S^ 


= [dS^/dfUf^ 


-1) 


so 






^J.- ^J.l 


= G-'''g-'S 










G{fA^ 


Gi + G[{fi-fi, 


)^ 


Gi 


[l + c 


iiSf 


where 












1 


G',G',"'g-' 











dfi \dfij Gi{l + aiSy 

We have derived this equation near fi = fii but away from there J^" is small and as dS^/dfi is small and 1 + aiS is not near 
zero, the equation can be taken to hold everywhere. We now integrate llAlOII and remembering that df/dfi is zero at fii we 
find 

^'^/'^^^ = G^lfa.Sr (^^^^ 

It was to get this simple result that we chose the cubic form in IIA8II . Taking the square root and multiplying up by 1/S + ai 
we find 

""f +«i(l-/) = Gr'^'(M-Mi) (A12) 

To proceed further we must integrate df/S{f). As the major item in (lAH is /'^""^ we shall assume that near the maximum 
of /, S^ can be approximated as Gi (1 — f^'^)- We do this by taking 

F= i[dln[S'(0)-S'']/dln/]j^i (A13) 
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and 

Cl = gGi/V (A14) 

We now use the method of paper I to give us j{l — f'^")^^'^df for large 77. 

Set / = 1_- ilng; d/ = -i^. Then J^" = (1 - ^L Inq^f' ^ exp-(ln52) = i/g^ so 

/(I - fT'^^'^df = i /(g^ - l)"^'^^rf'7 = =ch''^q so setting g = chA we find 

/ = 1 - = ln(c/iA) (A15) 

and from lA12f A is given by 

A + aiCi \n{chA) = VCiG^^^'\fi - m) (A16) 

and we remember that ai = ^G'lG^ g"^ — |^i72ir~T- To complete the solution we must impose the boundary conditions 

and so find the value of C needed to make / one at its maximum. However the boundaries vary from one application to 
another so it is time to consider the differential equations and their ranges one by one. 

The application in this paper is to the solution of equation l|88|l Xd^f/dX^ — —C'^vf^'^~^. 

To use the general theory above for this equation we write A for /i and G{X) — X, G\ = Ai and 5^ = C^(l — J^") the 
boundaries are at A = and A = 1 and V = v,G = Ci. At those boundaries we need chA = e" so that / = 0. If Ai is the positive 
value of A satisfying this then the negative solution is — Ai so at A = 1 and A = respectively Ai + aiCv — vCX^ (1 — Ai) 
— Ai + a\Cv = —vCX^ where ai — — r-m ■ 

Subtracting 2Ai = vCX~^''^ 

adding a\Gv = vGX^ (\ — Xi) = 2Ai(i — Ai) so using the value of ai above and eliminating C 

Ai(^-Ai)=J^/(12A?) 

Now Ai ~ !/ or more accurately i^ I 1 H J J so the quantity on the right is {-^) which is small so Ai is nearly 

a half. Accurately 



--\ 



1 + ^/1-^ 



(A17) 



with Ai known 

C = 2Ai''^Ai/i. (A18) 

This completes the solution when v is large 

/ ~ 1 - i In jc/i {2Ai(A - Ai) - ^{lnc/i[2Ai (A - Ai)]}| | (A19) 

AI Solutions virhen v is no longer large 

When V is lowered the solution will still have a single hump and will still fall to zero at the end points but it will no longer 
have precisely the form given in lA15II . IIA16t because there will be more deceleration of / away from /i = pi where our 
approximations for G(^) and S{f) are no longer valid. To allow for this we generalise the form of 1 A 1611 and <A17l l to 

f =1- ^\uchA (A20) 

A + A[nchA = hCG^^^^{fi-fii) (A21) 

where H and h can now be different and only for v large will they become equal to v. Furthermore A is no longer restricted 
to the form it took previously in terms of ai, i/ and C. We show in this section that with the four parameters H, A, h G and 
fii we can ensure that the solution has the right position of the maximum, the right deceleration there and that the gradients 
of the solution at the end points are appropriate for the differential equation. Since the desired solution has just one hump 
and is zero at the end points it is no surprise that we can fit it well not only for u large but for i^ > 1. In practice we find 
errors from computed solutions are less than 1% rising above 2% for u = 1. Since C is itself to be determined from the fitting 
of boundary conditions in the original equation we have actually five constants H, A, h, G and fii to be determined so we 
impose the following five conditions on our proposed form of solution A 

1/ & 2/ / = at the end points so there chAi = e" (A22) 
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(A23) 



„, (ff 1 fdS^\ ^,. . h^C^ 1 dS^ 

'/ - d^ = - 2g7 [wJ ,^, '^" """ ^- = 2 IT 

As we wish our proposed solution to solve IIA4fl as nearly as possible we now integrate ljA4^ by parts first from ^i to 1 and 
then from the lower boundary fim to /ii. These give us 

G(l) (^) - C (t-^ G\t.)dt, = [S{Q)f (A24) 



dii j j \ dn 

and 

"^^^'"^(i) +P (^) G'{p)df,^[Sm' (A25) 

We now use the forms 1A20II and <A21l l and substitute 

dfV 1 ,2,dA hC , , , 

J\ = th^A A26 

dnj H dii g\^^[1 + AthK\ 

into the integrals in 1A24L1A25II but use for the end value in 1A25II 

For (/')^^i we merely reverse the sign of A in ljA27^ . We then perform the integrals, so ljA22^ - ljA25fl provide us with the five 
equations to determine our five parameters. 

We wrote ljA26|l in that form because in the application to narrow jets G is linear so that the integrals are readily 
performed using A as the variable of integration. For that case the integral involved is 

-dA= / — dA = - — (1 + At/iA) ' ' y ^ ' 



l + AthA J 1 + AthA A' 'J (l + A)e2A + i-A 

-1 , /, A , .^ A A , /1 + A „2a\ 

= -ln(l + At/.A) + ^-^-^-^ln(^-^ + e -) 

Applying these methods to the solution of equation l|88|l . with A written for fj. and G{X) — X, S^ — C^(l — J^") and 
boundaries at A = and A = 1 the five equations are with Ai = ch~^e^ 

Ai+Alnc/iAi =/iCA7^/^(l- Ai) (A28) 

-Ai+Alnc/iAi = -/iCA+^''^ (A29) 

h'^/H = u (A30) 

and from 1A24II and llA25l l 



^2_ 1,, -2H, h^C hCX^f 1 



Ai A , / l + \/r"~e 



i+Av/r3-iFj^_i^__^i,( ;YrMJ KA3i) 



In solving these one uses the sum and differences of 1291 and 1301 1 and thinks of H rather than i/ as the independent variable. 
Figure 5 demonstrates the accuracy of the results compared with computed solutions. 

Finally we rewrite the approximate results as power series in u~^ 
H ^u- 1/42 - 5/26. i/'^ 

A = 1/2.1/"^ - 1/9.1/"^ + 1/12m-^ 

Ai = 1/2 - 7/31.1/"^ + 4/41.I/-2 

C, = hCX~^^^ ^2iy + \n4- 4/15.1/-^ 



